#include "StdAfx.h"
#include "PatternGenerator.h"

namespace SpikingNeuronLib
{
	// -------------------------------------------------------------------------------------------
	// Generate a firing pattern in which a single neuron is randomly selected (with replacement) to fire in each timeslot
	// Only one neuron can fire in each timeslot, therefore the average firing frequency per neuron cannot be greater than 1 Hz

	SingleRandomBackgroundPatternGenerator::SingleRandomBackgroundPatternGenerator(int totalNeurons)
		: PatternGenerator(totalNeurons)
	{
		_neuronFiringData = gcnew array<array<int>^>(1000);	// a slot for each msec

		Reset();
	}

	void SingleRandomBackgroundPatternGenerator::Reset()
	{
		ClearNeuronDataArray();

		// randomly select a neuron to fire in each one msec timeslot
		for (int i = 0; i < 1000; i++)
		{
			array<int>^ result = gcnew array<int> { _random->Next(_totalNeurons) };
			_neuronFiringData[i] = result;
		}
	}

	array<int>^ SingleRandomBackgroundPatternGenerator::Select(int msec)
	{
		return _neuronFiringData[msec];
	}

	// -------------------------------------------------------------------------------------------
	// Generate a firing pattern that is constant across resets
	// 'firingTimes' determines the msec timeslots in which firing occurs
	// and hence determines the overall network-wide firing frequency
	// The neurons that fire in available timeslots are randomly selected with replacement (excluding neurons in the 'inputPattern')
	// Only one neuron can fire in each timeslot, therefore the average firing frequency per neuron cannot be greater than 1 Hz

	FixedBackgroundPatternGenerator::FixedBackgroundPatternGenerator(int totalNeurons, List<int>^ firingTimes, Pattern^ inputPattern)
		: PatternGenerator(totalNeurons)
	{
		_neuronFiringData = gcnew array<array<int>^>(1000);	// a slot for each msec

		this->_inputPatternNeurons = this->CreateInputNeuronSet(inputPattern);		// may be set to nullptr

		this->SetRandomFiringTimes(firingTimes);

		Reset();
	}

	void FixedBackgroundPatternGenerator::Reset()
	{
		ClearNeuronDataArray();

		// randomly select a neuron to fire in each one msec timeslot
		for (int i = 0; i < 1000; i++)
		{
			array<int>^ result = nullptr;

			if (_randomFiringTimes == nullptr || _randomFiringTimes->Contains(i))
			{
				int neuron;

				if (this->_inputPatternNeurons)
				{
					// randomly fire only those neurons that are NOT in the input pattern
					do
					{
						neuron = _random->Next(_totalNeurons);
					}
					while (this->_inputPatternNeurons->Contains(neuron));
				}
				else
				{
					neuron = _random->Next(_totalNeurons);
				}

				result = gcnew array<int> { neuron };
			}

			_neuronFiringData[i] = result;
		}
	}

	array<int>^ FixedBackgroundPatternGenerator::Select(int msec)
	{
		return _neuronFiringData[msec];
	}

	void FixedBackgroundPatternGenerator::SetRandomFiringTimes(List<int>^ firingTimes)
	{
		if (_randomFiringTimes == nullptr)
		{
			_randomFiringTimes = gcnew HashSet<int>();
		}

		_randomFiringTimes->Clear();
		_randomFiringTimes->UnionWith(firingTimes);
	}

	HashSet<int>^ FixedBackgroundPatternGenerator::CreateInputNeuronSet(Pattern^ inputPattern)
	{
		HashSet<int>^ inputNeuronSet = nullptr;

		if (inputPattern)
		{
			List<FiringEvent>^ patternList = inputPattern->EventList;

			// create a set of all pattern neurons
			inputNeuronSet = gcnew HashSet<int>();
			for (int i = 0; i < patternList->Count; i++)
			{
				int neuron = patternList[i].NeuronIndex;
				inputNeuronSet->Add(neuron);
			}
		}

		return inputNeuronSet;
	}

	// -------------------------------------------------------------------------------------------
	// Generate a firing pattern from a Poisson process
	// The average firing frequency per neuron can be specified and the firing pattern changes across resets
	// The Poisson process determines the msec timeslots in which firing occurs
	// The neurons that fire in available timeslots are randomly selected with replacement (excluding neurons in the 'excludedNeurons')
	// Multiple neurons can fire in each timeslot

	PoissonProcessPatternGenerator::PoissonProcessPatternGenerator(int totalNeurons, int frequency)
		: PatternGenerator(totalNeurons), _frequency(frequency)
	{
		this->Initialise(totalNeurons, frequency, nullptr);
	}

	PoissonProcessPatternGenerator::PoissonProcessPatternGenerator(int totalNeurons, int frequency, List<int>^ excludedNeurons)
		: PatternGenerator(totalNeurons), _frequency(frequency)
	{
		this->Initialise(totalNeurons, frequency, excludedNeurons);
	}

	void PoissonProcessPatternGenerator::Reset()
	{
		this->CreateOneSecondRaster(_offsets, _frequency);
	}

	array<int>^ PoissonProcessPatternGenerator::Select(int msec)
	{
		array<int>^ result = nullptr;

		if (_neuronFiringData[msec])
		{
			result = _neuronFiringData[msec];
		}

		return result;
	}

	void PoissonProcessPatternGenerator::Initialise(int totalNeurons, int frequency, List<int>^ excludedNeurons)
	{
		_neuronFiringData = gcnew array<array<int>^>(1000);	// a slot for each msec

		if (excludedNeurons != nullptr)
		{
			_excludeSet = gcnew HashSet<int>();
			_excludeSet->UnionWith(excludedNeurons);
		}

		_offsets = gcnew array<int>(totalNeurons);
		for (int i = 0; i < totalNeurons; i++)
		{
			_offsets[i] = 1000;
		}

		this->CreateOneSecondRaster(_offsets, frequency);
	}

    inline int PoissonProcessPatternGenerator::NextISI(int freq)
	{
        return int (Math::Round(Math::Log(_random->NextDouble()) / (-(float)freq / 1000.0)));
	}

    // create a sequence of spike times separated by random ISI values
    void PoissonProcessPatternGenerator::CreateOneSecondRaster(array<int>^ offsets, int freq)
	{
        Dictionary<int, List<int>^>^ spikeDictionary = gcnew Dictionary<int, List<int>^>();    // spiking neurons keyed by time
		if (freq > 0)
		{
			for (int neuronIndex = 0; neuronIndex < offsets->Length; neuronIndex++)
			{
				if ((_excludeSet == nullptr) || (!_excludeSet->Contains(neuronIndex)))
				{
					int offset = offsets[neuronIndex] - 1000;
					// generate row
					while (offset < 1000)
					{
						offset += this->NextISI(freq);
						if (offset >= 0 && offset < 1000)
						{
							offsets[neuronIndex] = offset;    // last valid offset

							List<int>^ currentList;

							if (spikeDictionary->ContainsKey(offset))
							{
								currentList = spikeDictionary[offset];
							}
							else
							{
								currentList = gcnew List<int>();
								spikeDictionary->Add(offset, currentList);
							}
							currentList->Add(neuronIndex);
						}
					}
				}
			}
		}

		ClearNeuronDataArray();

		if (spikeDictionary->Count > 0)
		{
			// transfer the firing data to an array of arrays
			for each (KeyValuePair<int, List<int>^> kvp in spikeDictionary)
			{
				int index = kvp.Key;

				List<int>^ firedNeuronsList = kvp.Value;
				if (firedNeuronsList->Count > 0)
				{
					_neuronFiringData[index] = firedNeuronsList->ToArray();
				}
			}
		}
	}
}
